Limits of sensing temporal concentration changes by single cells 
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Berg and Purcell [Biophys. J. 20, 193 (1977)] calculated how the accuracy of concentration sensing 
by single-celled organisms is limited by noise from the small number of counted molecules. Here we 
generalize their results to the sensing of concentration ramps, which is often the biologically relevant 
situation {e.g. during bacterial chemotaxis). We calculate lower bounds on the uncertainty of ramp 
sensing by three measurement devices: a single receptor, an absorbing sphere, and a monitoring 
sphere. We contrast two strategies, simple linear regression of the input signal versus maximum 
likelihood estimation, and show that the latter can be twice as accurate as the former. Finally, 
we consider biological implementations of these two strategies, and identify possible signatures that 
maximum likelihood estimation is implemented by real biological systems. 



Cells are able to sense concentration gradients with 
high accuracy. Large eukaryotic cells such as the amoeba 
Dictyostelium discoideum and the budding yeast Saccha- 
romyces cerevisiae can sense very shallow spatial gradi- 
ents by comparing concentrations across their lengths . 
By contrast, small motile bacteria such as Escherichia 
coli detect spatial gradients indirectly by measuring con- 
centration ramps (temporal concentration changes) as 
they swim |2J, and can respond to concentrations as low 
as 3.2 nM — about three molecules per cell volume [3]- 
The noise arising from the small number of detected 
molecules sets a fundamental physical limit on the ac- 
curacy of concentration sensing, as originally shown in 
the seminal work of Berg and Purcell [H [S]. This ap- 
proach was recently extended to derive a fundamental 
bound on the accuracy of direct spatial gradient sensing 
[6|. However, no theory exists for the physical limit of 
ramp sensing, which is what bacteria actually do when 
they chemotact. In this Letter, we present such a theory 
for different measurement devices, from a single recep- 
tor to an entire cell. We contrast two strategies: linear 
regression (LR) of the input signal (in line with Berg 
and Purcell) and maximum likelihood estimation (MLE) 
[21 [8], a method from statistics to optimally fit a model to 
data, revealing an up to twofold advantage for the latter. 
Finally, we introduce a biochemical signaling network, 
similar to the E. coli chemotaxis system, that outputs an 
estimate of the ramp rate. Consistent with the derived 
theoretical bounds, we find that a mechanism emulating 
MLE yields twofold higher accuracy that one emulating 
LR. However, this improved performance has a cost: ei- 
ther storage of signaling proteins near the receptors, or 
irreversibility of the receptor cycle with concomitant en- 
ergy consumption. 

Sensing small numbers of molecules implies rela- 
tive noise ^n^^^^, where n is the number of detected 
molecules. Berg and Purcell (BP) calculated how this 
noise affects the accuracy of concentration sensing [4]. 
They considered three types of measurement devices: a 
single receptor, a perfectly absorbing sphere, and a per- 




FIG. 1: Schematic of measurement devices and correspond- 
ing time traces for linearly increasing concentration c{t) = 
Co + cit. (a) Left: a single receptor binds a particle at rate 
k+c{t), and releases it at rate k-. Right: binary time series of 
receptor occupancy, (b) Left: particles are incident on an ab- 
sorbing sphere with average flux 4:TrDac{t). Right: sequence 
of times when a particle hits the sphere, (c) Left: a moni- 
toring sphere counts the number of particles inside its volume 
without hindering their diffusion. Right: number N{t) of par- 
ticles inside the sphere as a function of time. 



fectly monitoring sphere. Following their approach, we 
investigate ramp sensing by these three devices when pre- 
sented with a concentration c{t) = cq + cit, as schema- 
tized in Fig. [l] 

A single receptor [Fig. [ija)] binds particles at rate 
fc_|_c(t) and unbinds them at rate fc_. Following BP, we 
assume that diffusion is fast enough that the receptor 
never rebinds the same particle. An ideal observer has 
access to the binary time series s{t) of receptor occu- 
pancy between —T/2 and T/2. The lengths of bound and 
unbound invervals have exponential distributions with 
means l/fc_ and l/Zc+c, respectively. Throughout, we 
assume that the ramp is shallow, ciT <^ cq, and that the 
observation time is long compared to receptor kinetics, 
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T ^ l/A:+c. In BP, the true concentration c is es- 

timated from the fraction of time the receptor is bound, 
s — ^ which is equal to the equilibrium 

occupancy in the limit of large times: 

s w (s) fc+c/(fc_ + fc+c), (1) 

where (•) represents an ensemble average. Following a 
similar strategy, we can estimate the ramp rate by per- 
forming the linear regression of s{t) to sq + sit: 
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dtts{t), (2) 



from which the concentration and the ramp rate are es- 
timated using (A18) as: 
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The uncertainties of these estimates can be calculated 
from the time correlations of receptor occupancy (see Ap- 
pendix A 1 a I , yielding: 

eg n' (co/r)2 
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where n is the total number of binding events in the time 
T. Note that the result for cq is precisely that of BP 

HIH]. 

In [S] , it was shown that the accuracy of concentration 
sensing could be improved using maximum likelihood es- 
timation. In this scheme, the parameters of the model are 
chosen to maximize the probability ("likelihood") that 
the observed data was generated by the model. Can we 
also improve the accuracy of ramp sensing over LR by 
using this method? The time trace s{t) can be char- 
acterized by the series of binding {tf ) and subsequent 
unbinding {t~) times, i = The probability of 

the data within our model is [8j: 



P : 



k^_llk+c{tt), (5) 

1=1 

where Ti, is the total bound time. The concentration 
and the ramp rate, Cq and ci, are the model parameters. 
Given the times of the events, the likelihood is maxi- 
mized with respect to cq and ci by solving dP/dco = 
and dP/dci — 0, from which the maximum likelihood es- 
timate {c^^^,cY^^) is obtained. In general these equa- 
tions have no simple solution, but we can obtain the aver- 
age behavior by exploiting the fact that binding and un- 
binding are fast with respect to concentration changes, 
i.e. that the receptor remains adiabatically in equilib- 
rium with the concentration c{t). We can thus simplify 
the sum and product in ([5|: 
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dtc{t) 



logc(t^ 



dt[l-(s(t))]c(i), (6) 



dtk_{s{t)) log c{t), (7) 



where (s(t)) is the equilibrium occupancy at time t, 
given by (AI8) with c = cq + cit, where cq and ci are 



the true parameters that generated the data. Applying 
this approximation to dP/dco, dP/dci, we confirm that 
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ci for T — T' 00 (see Appendix 
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AlbL For finite times, the errors in c^^^,cf^^ 



be estimated by the Cramer-Rao bound [5] , which states 
that the variance of parameter estimates exceeds the in- 
verse of the Fisher information, and approaches equality 
in the limit of long time series: 



{dc^dc)>^[dJd,iogP] 



(8) 



where <5c = (c^^e _ gg^cfLE _ 5^) and = 
{d/dcQ,d/dci). Again we can use the adiabatic approx- 
imation to compute the Hessian of the log-likelihood on 
the right-hand side of (A21), to obtain: 
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These variances are half the ones obtained from LR Q. 
The first result for constant concentrations is that of [8] . 
As observed there, the LR estimate adds the uncertain- 
ties from both bound and unbound interval durations. In 
contrast, the maximimum likelihood estimate relies only 
on unbound interval durations, since these carry all the 
information about the concentration. 

We now turn to ramp sensing by an entire cell, starting 
with the case of an idealized absorbing sphere [Fig.[ljb)]. 
An ideal observer witnesses a time series of absorption 
events, described by the instantaneous current /(<) = 
X]"=i ^(^ — tj), where 5{t) is the Dirac delta function 
and {ti} are the absorption times. The average cur- 
rent of molecules impinging on the sphere is given by 
{I{t)) — AnDac{t), where D is the diffusivity, a the 
sphere radius and c(t) the concentration far from the 
sphere [ij. Applying the same methods used for the 
single receptor, we calculated the uncertainty of ramp 
sensing for lin ear r egression of I{t) as well as for MLE 
(see Appendix A 2 1 . We found no difference between the 



two strategies, which both yield the same uncertainties 
as m with n now the total number of molecules ab- 
sorbed during time T: n ^ AirDacoT. For a monitoring 
sphere [Fig. [ijc)], molecules are free to diffuse into and 
out of the sphere, and the observer records the number 
N{t) of particles inside the sphere as a function of time. 
On average this number is {N{t)) — (4/3)7ra'^c(t). Per- 
forming a linear regression of N{t) to Nq + Nit, one can 
estimate the concentration and the ramp rate through 
c^" :== 3iVo/47ra3 and c^'" := mi/ina^. Following g], 
the uncertainty of these estimates can be calculated from 
the time autocorrelation of N{t) (see Appendix |A 3[ ), 
yielding: 
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The first result was obtained in g] . Maximum likelihood 
is difficult to implement in the context of the monitoring 
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FIG. 2: Biochemical network for measuring concentration 
ramps. Binding of ligand to the receptor increases its activ- 
ity u and causes species x to be produced. This production is 
downregulated by a feedback factor y which is itself catalyzed 
by X. Right: average network response to a step function in 
the concentration, c{t) = co + Ac9{t — to) (solid curves) and 
to a ramp, c{t) = co + ci{t — to)0{t — to) (dotted curves). In re- 
sponse to the step function, the network adapts precisely and 
X decays back to its original value after an initial increase. In 
response to a ramp, x shifts by an amount proportional to the 
ramp rate. The quantitative ability of the network to sense 
such ramps depends on whether receptors signal continuously 
or in a discrete burst upon ligand binding. 



sphere because it requires a sum over all possible histories 
of particles exiting and returning to the sphere. Thus, 
whether the LR result can be improved upon remains an 
open question. 

Maximum likelihood estimation is in general the op- 
timal way to sense ramps, and provides a twofold im- 
provement over simple linear regression in the case of 
the single receptor. Could MLE be implemented in bio- 
logical systems? To address this question, we now intro- 
duce a simple, deterministic biochemical network (Fig.[2| 
that can approach the optimal performance limit set by 
MLE. The same network implements either LR or MLE 
depending on the receptor signaling mechanism: LR is 
implemented if each receptor signals continuously while 
a particle is bound; MLE is implemented if each receptor 
signals with a fixed-size burst upon binding a particle, 
and then releases the particle rapidly. The first case cor- 
responds to integrating the fraction of time the receptor 
is bound, while the second corresponds to counting bind- 
ing events. Accordingly, we will show that the shot noise 
(Poisson noise) due to the stochastic nature of binding 
and unbinding is twice as large in the first case as in 
the second. Let u{t) be the receptor activity, propor- 
tional to the instantaneous production rate of signaling 
molecules. For continuous signaling, this activity is sim- 
ply proportional to receptor occupancy: u(t) = as(t), 
whereas for burst signaling, u(t) is a series of fixed-size 
bursts at the times of binding: u{t) — f3^" S{t — tf). 
Without loss of generality, we set a = fc_ and /3 = 1 so 
that {u{t)) is equal to the mean rate of binding events in 
both cases, {u{t)) = k-k+c{t)/{k_ + k^c{t)). For av- 
eraging times much longer than l/fc_ and 1/fc+c, we 



can approximate the fluctuations of u{t) by Gaussian 
white noise, u{t) = {u{t)) + Su{t), where {Su{t)du{t')) = 
g{u{t))S{t - t')/[l + k+c{t)/k_Y, with 5 = 2 for con- 
tinuous signaling, and 5 = 14- [k^cjk^Y for fixed-size 
burst signaling (see Appendix B 1 ) . For rapid unbinding, 
fc_ — >■ -|-oo, we recover the same twofold difference as be- 
tween ^ and (|9]), and for the same reason: in the case 
of continuous signaling, noise from the stochasticity of 
bound intervals adds to the noise from random arrivals. 

To extract the ramp rate from receptor activity re- 
quires a network that "takes the derivative" of its input 
signal. An example is the E. coli chemotaxis system, 
which relies on precise adaptation via integral feedback 
[TOl [11] . A minimal deterministic version of such a net- 
work is schematized in Fig. [2] and described by the fol- 
lowing differential equations: 



dx 
~dt 



kx [uf{y) - x] 



dy 
dt 



kyiX-l); (11) 



where for simplicity u{t) is the activity of a single re- 
ceptor, and X is the concentration of signaling molecules 
it produces, /(y) is a monotonically decreasing func- 
tion regulating the production of x. The role of y is 
similar to that of the receptor methylation level in E. 
coli: y precisely adapts the production rate of signal- 
ing molecules so that the steady-state value of x does 
not depend on the external ligand concentration. This 
property is illustrated by the graphs on the right side of 
Fig. [2] which show how the network responds to a sud- 
den change in ligand concentration (solid curves). While 
the network output x is insensitive to the absolute con- 
centration, it responds to steady ramps (dotted curves). 
When the input varies slowly in time, {u{t)) — uq -\- uit 
(with ui <ti uokxjUoky), the system responds by shift- 
ing X away from 1 so that the change in y{t) tracks the 
change in u(t): 

(^W> = i + 7t^, {y{t))^yo-j'^~j-t, (12) 

kyUo kyUo uq 

with uofiyo) = 1 and 7 = - f (yo) / f (yo) ■ Thus y pro- 
vides a readout of the absolute concentration, and x pro- 
vides a readout of the ramp rate. The accuracy of these 
representations is limited by the ligand binding shot noise 
6u{t). The effect of noise can be calculated by expand- 
ing the solution of (111 linearly around its average (see 
Appendix B 2 ) : 



" Sx{t) ' 




_ Sy{t) 





x{t) - {x{t)) 

y{t) - {y{t)) 



dt'K{t~t')Su{t') 



with K(<) = — e 



Uq 



cosh(wt) - |j sinh(wt) 
— sinh{ujt) 



where — — k^ky/j {uj can be imaginary). From 
( 12 1 we deduce the uncertainties of co and ci: 
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For a fixed ky, the optimal value of k^ is the smallest one 
with a non-oscillating response kernel K(t): k^ = 4fcy/7. 
Systems with oscillating kernels are undesirable because 
they detect oscillations rather than ramps. For k^ = 
Aky/^, our results are consistent with those of Eqs. Q 
and (|9]), namely uncertainties inversely proportional to 
the number of binding events, if we interpret j/ky — > T 
as the effective time of measurement, and uq as the rate 
of binding events. The factor g reflects the difference 
between the two mechanisms of receptor signaling. 

Despite its simplicity, our biochemical model may help 
analyze features of real biological systems. There are two 
separate aspects to the model: on the input side, different 
mechanisms of receptor signaling — continuous signaling 
(LR) versus burst signaling (MLE) — affect readout ac- 
curacy; on the output side, integral feedback provides a 
natural readout for sensing ramps. 

Many receptors, including the well-studied chemotaxis 
receptors of E. coli, signal continuously rather than in 
bursts, and therefore do not employ MLE. In practice, 
how could cells implement MLE? A receptor could simply 
"store" a fixed amount of signaling molecules and release 
all of them upon ligand binding. Alternatively, receptors 
could signal continuously following a binding event but 
with a narrowly peaked distribution of durations. Our re- 
sults can easily be extended to an arbitrary distribution 
of durations yielding g = 1 + {{Sti,)'^) / {rt)^ : the more 
peaked the distribution of Tf,, the less noisy the readout. 
For equilibrium binding/unbinding, we find g > 2 (see 
Appendix [C| , with an irreversible binding cycle driven 
by energy dissipation required to achieve g < 2. Inter- 
estingly, there are examples of such irreversible cycles in 
ligand-gated ion channels [T2HT4] . where ions play the 
role of our output signal x. In these ion channels, peaked 
open-time distributions are interpreted as evidence that 
time reversibility is broken and energy is being consumed 
[l5] . We speculate that the role of this irreversibility may 
be to reduce the variance of bursts, thereby increasing the 
accuracy of concentration or ramp sensing. Relatedly, a 
multiplicity of irreversible steps in rhodopsin signaling 
has been shown to explain the reproducibility of single- 
photon responses in rod cells [TB]. 

As for the mechanism of ramp sensing, the integral 
feedback system underlying E. coli chemotaxis is similar 
to our simple model. However, the receptor mcthylation 
level, which plays the same role as y in our model, adjusts 
the binding/unbinding rates k+/k- so that fc_ ~ fc+c, 
rather than adjusting the production rate kxf{y)u as 
in (111. h\ E. coli receptors increase their gain by re- 
sponding cooperatively (17j . and k- ^ k+c is required to 
maximize this gain, which precludes the limit k^ k^c 
required for MLE. Moreover, k^ is physically limited by 
diffusion and receptor size, and should optimally be kept 
near the diffusion limit to maximize the number of bind- 
ing events. It is worth noting that in E. coli the mcthy- 
lation and demethylation processes responsible for inte- 
gral feedback are themselves subject to noise, giving rise 
to additional fluctuations [TB]. For a receptor signaling 



in bursts, integral feedback could act by adjusting the 
number of released molecules upon binding if the recep- 
tor stores molecules, or the mean bound duration (rf,) 
if signaling is continous, or the channel conductivity in 
ligand-gated ion channels. We hope that our analysis will 
suggest experiments for testing these scenarios. 

We thank Pankaj Mehta and Aleksandra Walczak for 
helpful suggestions. T. M. was supported by the Human 
Frontier Science Program and N.S.W. by National Insti- 
tutes of Health Grant No. ROl GM082938. 



Appendix A: Uncertainties in ramp sensing 

1. Single receptor 

We consider a single receptor, which particles bind to 
and unbind from stochastically. The unbinding rate is 
denoted by fc_ , and the binding rate by k+c{t), where c(t) 
is the concentration of particles. We assume fast diffusion 
and so neglect rebinding of particles f5]. Between the 
times — T/2 and T/2, the concentration changes slowly 
with time: c{t) = cp + Cit. We assume that variations are 
small: ciT <C cq, and that the integration time T is large 
compared to the waiting times: T 3> fcZ^, (fc+co)^^. 



a. Linear regression 

We perform linear regression of the binary time series 
s{t) of receptor occupancy s{t) to Sq + Sit. The time 
trace s{t) can be characterized by the series of binding 
{tf ) and subsequent unbinding times, i = 1, . . . ,n. 
The resulting estimates for co and ci are: 
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^ j^^^^dtts{t)^^Y.^tT^i ),(A4) 



are the results of the linear regression. 

In the limit of long time series, these estimates give the 
correct answers: replacing s{t) by its ensemble average 
value, (s(i)) = fc+(co + ci<)/[fc_ + fc+(co -f cit)], where 
co,ci are the true values of the concentration and the 
ramp rate, one obtains: 



(so) 
(si) 



-Co 



k. 



k+co ' 
-k+ci 



(fc_ + k+coy 
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which leads to Cg^ = cq and = Ci on average. The 
expected error can be obtained from the covariance ma- 
trix (Js*(5s). To compute this covariance matrix we need 
the following quantity: 

(.W,«+t'))-(.(0>(.(t+i'» - - V +^,c(0]^ ■ 

(A7) 

where the approximation « is valid provided that ^ 
cq/ci. In the limit of large times we have: 



{SsqSsi) 



2 fc_fc+co 

r(fc_ + fc+co)3' 

24 k^k+co 
r3 (fc_ + A:+co)3 ' 
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from which we deduce (using Eqs. Al A2 and ciT ^ cq) 
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The above result for cq was originally obtained by Berg 
and Purcell g]. 



b. Maximum likelihood 



The probability ( "likelihood" ) of a sequence of binding 
and unbinding events at times {tf} and {t~} between 
-T/2 and T/2, is: 



P = 



e — "e 



{k^k+ri[citt). 



(A13) 

The log-likelihood is (up to a constant independent of cq 
and Ci): 

log P = const - ^+ E (^0 + c /'^'^^^' ) - 

+ Elog(co + cii+). 

(A14) 

Given the times of the binding and unbinding events, the 
optimal strategy for estimating the concentration and the 
ramp rate is to maximize this log-likelihood with respect 
to Co and ci, i.e. to solve 

^=0 and ^=0, (A15) 

OCo OCi 

from which we can obtain the maximum likelihood esti- 



mates Cq and Cj To solve Eq. A15 we exploit the 



self- averaging property of log P as T — -|-oo (adiabatic 
approximation). The sums in Eq. A13 and A14 become: 



n 



dt[l- {s{t))]c{t), (A16) 
dtk_{s{t)) log c{t),{A17) 



where {s{t)) is the equilibrium probability of the receptor 
being bound at time t: 



{s{t)) 



k+c(t) 
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and c(t) = Cq + cii is the true concentration. We take the 
derivatives of the log-likelihood (Eq. A14 with replace- 
ments from Eqs. A16 and A17) with respect to cq and 
ci: 



91ogP 
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fc_ + /s+(co + Cit) 

k-k+{co + cit) 
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-T/2 [k- + k+ (co -I- cit)] (co + cit) ' 

'^/^ k k,f 

dt f-f+^ , ^ (A20) 

-T/2 fc_ + fc+(co + Cit) 

k_k+{^a + c^t)t 
-T/2 [k- + k+ (cq + Cit)] (co + Cit) ■ 



It is straightforward to confirm that these expressions 
become zero for cq = Cq and Ci = Ci. Thus, in the limit 
of long time series, the maximum likelihood estimates 
coincide with the true values of the concentration and 
the ramp rate. The expected error of these estimates is 
given by the Cramer-Rao inequality, which becomes an 
equality in the limit of large T: 



{Sc'^Sc)>-[dJd^\ogP] \ 



(A21) 



do,cf^^ - ci) and = 



where 6c = (cjf^^ 
{d/dco,d/dci). 

The second derivatives of the log-likelihood are, to 
leading order: 



\ogP 



logP 



fc_fc+(5o + git) 

_T/2 [fc- + fc+(co + Cit)](co + Cit)2 

T k^k+ / . „„x 

— , ^, , (A22) 

Co K- + K+Co 
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logP 
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dt ■ 



k-k+{c.Q -\- Cit)t^ 



T/2 [fc- + fc+(co Cit)](co + Cit)2 
'+ (A23) 



12co fc_ + fc-i-Co ' 
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dt ■ 



k^k+{cQ + cit)t 



T/2 [fc- + fc+(co + Cit)](co + Cit)2 

(A24) 
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Exploiting the diagonal structure of the Hessian, we ob- 
tain: 



(co/T)2 
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where n is the total number of binding events in the time 
series, n « k-k^c^T / [k- + fc+co). The above result for 
Co was first obtained in [6]. 



2. Absorbing sphere 

We now consider as a measurement device a perfectly 
absorbing sphere of radius a. Solving Laplace's equation 
with absorbing boundary conditions at the surface of the 
sphere yields the average flux of particles impinging on 
the sphere [4j: 



(/(t)) =A-KDac{t), 



(A27) 



where D is the diffusivity, a is the sphere radius, and 
c(i) = Co + cit is the (time-dependent) concentration far 
away from the sphere. When the ramp arises from motion 
of the sphere up a spatial gradient (with velocity v), one 
must be careful: in general the flux of particles is not 
uniform around the sphere, and will differ between the 
front and the rear of the sphere (windshield effect) Hj. 
However, this effect can be neglected when measuring the 
change of concentration with time if the particle turnover 
rate D jo} is much larger than rate of moving one body 
length w/a, i.e. if Dja ^ v. For bacteria the relevant 
numbers are D k, 10"^ fiw? /s, a ~ 1 /zm, and w « 30 
/im/s, so that D/a ~ 10'^ /Ltm/s. Therefore, from the 
perspective of a swimming bacterium, the concentration 
appears to be changing with time, essentially uniformly 
in space, which is the case we consider. 



a. Linear regression 

For a constant concentration, the simplest estimate for 
the concentration is Cq^ = n/AnDaT, where n is the total 
number of absorption events in time T. In the presence 
of a ramp, an estimate of both concentration and ramp 
rate can be obtained by performing the linear regression 
of 



(A28) 



to AnDa^CQ + cit), where the {ti} are the absorption 
times, yielding: 
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^ S{t - ti) - AnDaico + at) 
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Since the particle absorptions are independent events, n 



is a Poisson variable, {{Sn)'^) = 



AttDqco (as above. 



co,ci are the true values of the concentration and the 
ramp rate), and thus ((^Co)^)/cq = 1/n. Concerning the 
ramp rate, we have 



n \ i-T/2 

"^tA ^ / dtt{I{t)) = AnDaci 

; 1 / <^—T/2 
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and 




\ dtt'^ {I (t)) = At: DacQ 
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where we have used the fact that over a short time St, 
the number of absorbed particles is a Poisson variable, 
which entails: 




t+st 



Therefore: 



dt' t'5I{t') 



{{Sen') = 




12co 



Finally, in summary: 
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12 

n 



I{t))t^6t. (A33) 



(A34) 



(A35) 
(A36) 



where n is the total number of absorption events in time 
T. 



b. Maximum likelihood 

The observations by an absorbing sphere in the interval 
(— T/2,T/2) are summarized by the sequence of times 
when particles are absorbed, {ti}, i = 1, . . . ,n. The log- 
likelihood of this sequence reads: 



/T/2 n 
d<(/(t)),„,,, +^log(/(<,))c„,c, 
'T/2 t" ; 



47rL>aT' 



const — AttDucqT + ^ log(co + citi), 

(A37) 



4=1 
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where {I{t))co,ci is the expected average flux for concen- 
tration Co and ramp rate ci. Setting the derivatives to 
zero, dlogP/dco = 0, d\ogP/dci ~ 0, yields to leading 
order: 



MLE 
^0 



MLE 



AirDaT ' 



MLE 



HA 



(A38) 
(A39) 



where we have used C\T ^ cq. In the limit of long time 
series, the sums self-average and become: 



^tiW^^t,y=y" dtATrDa{cQ+cit)t 



2^3 

12' 



(A40) 



(A41) 



T/2 

dtAnDa{co + cit)i^ 

(A42) 

"l2' 



A38 and 



which yields Cq^^ = So and cf^^ = ci using Eqs. 

EMI 

To obtain the uncertainties, we calculate the second 
derivatives of the log-likelihood and use the Cramer-Rao 
bound. After a calculation similar to Eqs. A22||A24" we 
find: 



1 
n 
12 

n 



(A43) 
(A44) 



Therefore, maximum likelihood estimation is equiva- 
lent to linear regression in the limit of long time series. 
(One can similarly show that the same conclusion applies 
to the sensing of spatial gradients: the estimate from lin- 
ear regression derived in is equivalent to maximum 
likelihood estimation.) 



3. Monitoring sphere 

Linear regression 

Following Berg and Purcell [4], we now consider a 
perfect monitoring sphere of radius a, which lets par- 
ticles freely diffuse in and out of the sphere. This device 
records the number of particles N{t) inside the sphere at 
all times, with ensemble average {N{t)) — (4/3)7ra'^c(i). 
The particle turnover time, /D (where D is the diffu- 
sivity), is assumed to be small compared to the total time 



T of observation. Regressing N{t) to (4/3)7ra^(co -|- cit) 
yields an estimate for the concentration and the ramp 
rate: 

fT/2 

' dtN{t), (A45) 

-T/2 



Co — 



1 1 

vT 



LR 



1 12 



.T/2 



= / dtN{t)t, 



-T/2 



(A46) 



where V — (4/3)7ra'^ is the volume of the sphere. In the 
limit of long times, this estimate gives the true value 
of the concentration and the ramp rate: 



„LR\ 



{Ci'-'') = ci. The uncertainties are given by the vari- 
ances of and c^^. For example, the uncertainty of 
the concentration estimate is: 



,LR,\2\ 



1 1 

y2j^ 



T/2 



dt 



T/2 



T/2 



T/2 



dt' {SN{t)SN{t')), 



(A47) 

where SN{t) = N{t) - {N{t)). To compute 
{SN{t)SN{t')), we decompose N{t) as a sum over M in- 
dependent particles in a large volume U ^ V: N{t) = 
X]i=i ^i(J')j where Xi{t) is a binary variable for each par- 
ticle, whose value is 1 if the particle is in the sphere, and 
otherwise. We assume that M and U are large, such 
that M/U = Co + cit. The ramp can arise from a change 
in volume U, or from creation/annihilation of particles. 
For definiteness, we will assume that the number of parti- 
cles M stays fixed, while the confining volume U varies in 
time. Since the particles are independent of each other, 
we have {6xi{t)Sxj{t')) = 0, and therefore: 



1 1 

Y2J^ 



T/2 



T/2 



T/2 



T/2 



dt I dt'{5x{t)5x{t')), 
(A48) 

where x{t) the binary indicator variable of a single par- 
ticle, and 5x{t) = x{t) — {x{t)). In the limit of large time 
differences |i — <'| ^ D/o?' , the autocorrelation fun ction 
{5x{t)5x{t')) decays to 0. Since T > D/a^, Eq 
simplifies to: 



A48 



LR\2\ 



2M 



1 1 



dt{5x(Q)Sx[t)) , (A49) 



where we have used time translation invariance and 
time-reversal symmetry between t and t' . The quantity 
{x{Q)x(t)) is the probability that a particular particle was 
in the sphere at times and t, i.e. P{x{Q) — l)P{x{t) = 
l|a;(0) = 1) = (a;(0))uBp(t), with {x{Q)) = (Sq -f 5it)t;/M 
and 



MBP(0 :=P(a:(t) = l|.T(0) = l) 



1 

V 



df I df 



exp 



|r— r I 



(A50) 



ADt 

(47rDi)3/2 ' 

where the integrals are over the sphere of radius a. A 
calculation familiar from electrostatic done in |4j yields: 

dtuBp{t) = --pr. 



(A51) 
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Combining the above results and using {x{Q)){x{t)) = The calculation of the uncertainty of the ramp rate pro- 
{cV/M)^, we obtain the Berg and Purcell bound: ceeds very similarly: 

4 ~ 5nDacoT- ^^^^^ 



I.T/2 ,TI1 

/ dt \ 

J-T/2 J-T/2 



dt'tt' {Sx{t)Sx{t')) « 2M 



rT/2 

/ dtt 

J-T/2 



2 (cO 



iit)V 



M 



+00 



dt' UBp{t') 



r 



(A53) 



where we have used t' t. This approximation is justi- 
fied by the fact that for \t-t'\ > D/a^, {5x{t)5x{t')) 0. 
Finally we obtain: 



(co/T)2 



36 

bnDacoT' 



(A54) 



Appendix B: Uncertainties in a ramp-sensing 
biochemical network 



1. Input noise 

We first consider the most general case of the biochem- 
ical network shown in Fig. 2 of the main text. Ligands 
bind the receptor at a rate k+c{t). When a ligand binds 
to the receptor, it remains bound for a time t^, and (3 
signaling molecules are released. The distributions of 
both Tb and /3 can be arbitrary. For continuous signaling, 
particles are produced at a rate /3 while the receptor is 
bound, yielding /3 = arb, while for fixed-size burst sig- 
naling, /3 particles are realeased directly upon binding at 
time tf. Both scenarios can be described by the instan- 
taneous receptor activity: 

Continuous signaling u{t) = as{t), 
Fixed-size burst signaling u{t) = /3^^(5(i — tf), 

where, in the first equation, ,s(/:) = 1 when the recep- 
tor is bound and otherwise. W^icn there is only one 
bound state, Tb has an exponential distribution, with 
{Tb) = 1/k- and {{STb)^) = For averaging times 

much longer than the durations of bound and unbound 
intervals, u{t) can be approximated by a Gaussian vari- 
able. The mean of u{t) is: 



Continuous signaling {u{t)) = 
Fixed-size burst signaling {u{t)) = 



ak+c{t) 

l/{Tb}+k+c{t)' 

in) + i/k+c{t) ■ 



To estimate the variance of u{t) for fixed-size burst sig- 
naling, we calculate the variance of the number of binding 
events during a time At ^ (Tb), l/k+c{t) (but At <C T). 



The variance of the duration of one binding/unbinding 
cycle is: (.5(t6+t„)2) = {{5Tby) + l/[k+c{t)]'^ . If there are 
n binding/unbinding events during At, the relative vari- 
ance of n for a fixed duration At is equal to the relative 
variance of At for a fixed number of binding/unbinding 
events n: 

{{6nr) _ {{SAtr) _ 1 {i5Tbr} + l/[k+c{tr 

„2 (A+\2 „ r/_ \ , 1 _/.m2 ' \^^> 



{Aty n [{Tb) + l/k+c{t)Y 



hence: 




t+Ai 



dt'u{t')\ )=l3^{{Snf) 



= B^At 



{i5Tbr) + i/[k+c{t)f 



(B2) 



[{Tb) + l/k+cit)f 

Thus, for averaging times much longer than (rft) and 
l/fc+c(i), we may write: 

2. ^,,,i + {{Snr)[k+c{tf 



{Su{t)Su{t')) = I3^k+c{t) 



[1 + {Tb)k+c{t)r 



S{t~t'). 
(B3) 

Similarly, we can estimate the fluctuations of u{t) for 
continuous signaling: 

k+c{t) 



{5u{t)Su{t')) = a'^{T^) 



[1 + {Tb}k+C{t)]3 



S{t-t'). (B4) 



Setting a = 1/ {Tb) and /3 = 1 without loss of generality, 
we obtain for both cases: 



{Su{t)5u{t')) 



k+c(t) 



l + {Tb)k+c{t)' 
k+c{t) 



(B5) 

-Mt-t'), (B6) 



ni + (r,)fc+c(i)]3 

with g = 1 + {{STb)^) / (Tb)^ for continous signaling, and 
g = 1 + {{STbf) [k+c{t)Y for fixed -size burst signaling. 
Note that for an exponential distribution of bound inter- 
vals durations (((^Tf,)^) = — so that 5 = 2 for 
continuous signaling, and g = 1 + [k^c{t) / k-]"^ for fixed- 
size burst signaling. In the limit of short bound-time 
durations, Tb — )• 0, the general results become: 

{u{t)) = k+c{t), {5u{t)5u{tf)) = g{u{t))6{t - t'), (B7) 

with g = 1 -\- {{STb)^) / (Tb)^ for continuous signaling and 
g = 1 for fixed-size burst signaling. 
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2. Output noise 

Consider the network shown in Fig. 2 of the main text 
and described by the equations: 



dx 

~di 
dy 

dt 



= kx [uf{y) - x] 

= ky{x-l). 



(B8) 
(B9) 



When the network is presented with a slow concentration 
ramp cq + cit (such that cit <C cg), with the input u{t) 
given by Eq. |B7[ i.e. u{t) = uq + uit + 6u{t) with: 

uo = fc+co/(l + (Tf,)fc+co), (BIO) 

Ml = fc+Cl/(l + (T6)fc+Co)^ (Bll) 



the average network response is: 

Ul 



S{t-t'), (B12) 



(xit)) = 1 + 7 

{y{t)) = yo-i 
"o/(yo) = 1 



kyllQ 

2 Ul 



Ul 

kyUo Uq 



(B13) 

(B14) 
(B15) 



with 7 = —f{yo)/f'{yo)- We collect the fluctuations of 

X and y, 5x{t) = x{t) — {x{t)) and 5y{t) = y{t) — {y{t)), 
into a single vector: 



X 



Sx 
6y 



(B16) 



Linearizing Eqs. (B8) and (B9) yields: 



— — ^MX = — 

dt Uq 



Su{t) 




, with M = 



-ky 



(B17) 

Multiplying by e**^ on both sides and integrating, one 
obtains: 



X 



Uo 



dt' e 



Su{t) 




(B18) 



with 



tM 



cosh(wt) + ^ sinh(wt) 



: sinh(wt) 



sinh(wi) cosh(a;t) - ^ sinh(wi) 



(B19) 



cj^ = — kxky/j if kx > Aky/j, and the same expression with sinh sin, cosh — )■ cos, = kxky/y — 
otherwise. The fluctuations of x are given by: 



6xit) 



Uo J- 



cosh(a;t) smh{ujt) I ( cosh(cL>t') H ^ sinh(cL>t') I H — sinh(a;t) sinh(cjt') 



2uj 



cosha;(i — t') — sinha;(t — t'] 

^ ^ 2uj ^ ^ 



2uj 



5u{t'), 



6u{t'), 



and the average uncertainty: 



{(Sxy) 



k+Co{l + (rf,)fc+Co) J_ 



cosha;(t — t') — sinhu;(^ — t') 
2uj 



gkx 



2/0+00(1 + (rf,)fc+Co) 



, (B20) 



where the last equality is valid at steady state. Exactly 
the same result is obtained when k^ < Aky/y. Using 
Eqs. |Bll] and |B13l we derive the uncertainty of the ramp 
rate readout: 



{cokyl-fY 



gkx 

2uo ' 



(B21) 



In addition we can also evaluate the variance of y: 

Sy{t) = ^ /* dt' e-'=-(*-*')/2^ sinhw(t - t')5uit'), 

Uo J -00 ^ 

(B22) 

and therefore: 
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k+co{l + {n)k+ca) 



^^,g-/c(t-t')[sinhLj(f-t')]^ = 



2fc+Co(l + (Tf,)fc+Co) ' 



(B23) 



from which we obtain the uncertainty of the concentra- precisely, 



tion readout, using Eqs. BIO and B15 



((M 



2uo 



(B24) 



Appendix C: Distribution of intervals of signaling 
activity for a receptor at equilibrium 

We now prove that for any receptor in thermal equi- 
librium: 



{{dry 



> 1, 



(CI) 



where r is the duration of an interval of uninterrupted 
signaling activity. With the results from the previous 
section, this proves g > 2 for any continuously signal- 
ing receptor, unless free energy is consumed in the bind- 
ing/unbinding cycle. Previously and in the main text, the 
receptor was assumed to signal whenever the receptor was 
bound, and the receptor was inactive otherwise (t is then 
equal to the bound duration r;,), but our proof allows 
for a more general definition of signaling activity. We 
model the receptor by a Markovian system with TV dis- 
tinct states. M of these states are called active and form 
the subset ^ C {1, . . . , N}. The remaining states belong- 
ing to the complementary subset / = {1, . . . , are 
called inactive. 

Our goal is to calculate the distribution of intervals of 
uninterrupted activity when the system is in equilibrium, 
i.e. intervals where the system remains in A. In particu- 
lar, we will show that this distribution is a weighted sum 
of exponentials with positive weights: 



{{dry 



o v^-A/ 1 



- 1 = 



{j2k=i 



= 1 ^fcA^ 



-1>1, (C3) 



where the last step follows from Jensen's inequality, 
which becomes an equality if and only if M — 1. 

We now prove Eq. C2 Let Q = {qij} denote the ma- 
trix of rates from state j to state i, with qu = — Iji- 
Call \p) — {pi, . . . ,pn)'^ the vector of equilibrium prob- 
abilities, which satisfies Q\p) — \p). Because of detailed 
balance, we have qijPj = qjiPi- With proper ordering 
of indices, Q can be decomposed according to the active 
and inactive states: 



Q = 



]aa Qai 

^lA Qll 



(C4) 



Note that in [T^, Qai and Qia are square, diagonal ma- 
trices. Assuming that the system enters an active state 
at time t = 0, this state will he i A with probability 



[QAl\Pl)]i 

{^a\Qai\pi 



(C5) 



where Qai is the submatrix of rates from inactive to 
active states, \pi) is the projection of \p) onto subset /, 
and 1a is a vector of ones of dimension M. Starting from 
an active state j at time 0, the probability of still being 
active in state i € A at time t (if i e A), or alternatively 
to have exited the active states via the inactive state i € I 
at time t' < t {ii i ^ I), is given by Xij, where X is the 
solution to the following Master equation: 



dX 
~dt 



= Q'^^'^X, 



X = e*«° 



(C6) 



M 



Wk > and 



fe=i 



M 

k=l 



(C2) 

This result generalizes the one proven in the support- 
ing information of [19 . In that proof, active states were 
paired one-to-one to inactive states (whence TV = 2M), 
and, among the possible transitions from active to inac- 
tive states, only the ones within these pairs were allowed. 
Our proof places no restrictions on the number of states 
or on transitions between states. 



A consequence of Eq. C2 is that the relative variance 
of intervals is bounded below by one. This lower bound is 
only attained for a single-step process (M = 1) for which 
the distribution of intervals is a pure exponential. More 



where Q is the same as Q but with absorbing inactive 



Wk^l. states: Qff = and QJf = 0. Thus 



Q 



abs 



Qaa 
Qia 



and 



QlAQAAi^ 



tQA/ 





1) 1 



(C7) 
(C8) 



Therefore the total probability of having exited the active 
states at time t' < t, after starting in j G A a,t t = 0, is: 



[{lllQlAQA^Aie'"^' 



1) 



(C9) 



where (1/| is a vector of ones of dimension N — M. Fi- 
nally, the total probability of having exited the active 
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state at time t' < t after entering it at time t = is: left eigenvectors of Qaa- P{t) can now be rewritten as: 

{1a\Qai\pi) ^^^> {Ia\Qai\pi) ■ ^ ' 

which yields the probability distribution of intervals of have' 
activity: 

p., dC.. {h\QjAe-^--QAi\Pi) (v'^'^IQaiIpi) ^ E ^Pj'^'im 

dt {Ia\Qai\pi) _ 2 

= 2^ "jPj qijPj 

To exploit the property of detailed balance, we sym- 

metrize Qaa by defining Q = {qijy/pj /Pi\teA,jeA- Q ^ 



E 



being symmetric, it can diagonalized in an orthonormal 
base: ^e/je^ 

Q = -^\i,\u^){u\ withAfc>0 and (u'=|u'=') = 4fe'- ■= ^ 

(C12) and therefore: 

The diagonal form of Qaa thus reads: 



Qaa = -5]A,>i^^'=)(,;L^'^1, (C13) EfcLi 



fc=i 

which establishes Eq. |C2| with zii^ 

„R.'= _ „,fc /;rTT _ „,fc / z;^ ^;„u+ /'/l2/\ \//Y^ 42 



where ' ' = ^'^'^ ""i ' ^ "i/VPi are the right and (^fc/Afc)/(I]fc' ^fe'/-^fc')- 



[1] R. A. Arkowitz, Trends Cell Biol 9, 20 (1999). 

[2] R. M. Macnab and D. E. Koshland, Proc Natl Acad Sci 

USA 69, 2509 (1972). 
[3] H. Mao, P. S. Cremer, and M. D. Manson, Proc Natl 

Acad Sci USA 100, 5449 (2003). 
[4] H. Berg and E. Purcell, Biophys J 20, 193 (1977). 
[5] W. Bialek and S. Setayeshgar, Proc Natl Acad Sci USA 

102, 10040 (2005). 
[6] R. G. Endres and N. S. Wingreen, Proc Natl Acad Sci 

USA 105, 15749 (2008). 
[7] A. Edwards, Likelihood (Cambridge University Press, 

Cambridge, 1972). 
[8] R. G. Endres and N. S. Wingreen, Phys Rev Lett 103, 

158101 (2009). 

[9] S. M. Kay, Fundamentals Of Statistical Signal Processing: 
Estimation Theory (Prentice Hall PTR, Upper Saddle 
River, NJ, 2001). 
[10] N. Barkai and S. Leibler, Nature 387, 913 (1997). 
[11] T. M. Yi, Y. Huang, M. I. Simon, and J. Doyle, Proc 



Natl Acad Sci USA 97, 4649 (2000). 
[12] R. Schneggenburger and P. Ascher, Neuron 18, 167 
(1997). 

[13] S.-Q. Wang, L.-S. Song, L. Xu, G. Meissner, E. G. 
Lakatta, E. Ri'os, M. D. Stern, and H. Cheng, Biophys J 
83, 242 (2002). 

[14] L. Csanady, P. Vergani, and D. C. Gadsby, Proc Natl 
Acad Sci USA (2009). 

[15] D. Colquhoun and A. Hawkes, in Single- Channel Record- 
ing, edited by B. S. B and E. Neher (Plenum, New York, 
1995), pp. 397-482. 

[16] T. Doan, A. Mendez, P. B. Detwiler, J. Chen, and 
F. Rieke, Science 313, 530 (2006). 

[17] R. G. Endres and N. S. Wingreen, Proc Natl Acad Sci 
USA 103, 13040 (2006). 

[18] T. Emonet and P. Cluzel, Proc Natl Acad Sci USA 105, 
3304 (2008). 

[19] Y. Tu, Proc Natl Acad Sci USA 105, 11737 (2008). 



